function trace_average()
%try
sf = getmainselection;
if(sf > 0) 
   pflag = getplotflag;
   QueMessage('Trace Averaging', 1); % clear the que
   hf = findobj('Tag', 'TRAVG');
   if(isempty(hf))
      hf = figure('Tag', 'TRAVG', 'Name', 'Averages Plots', 'NumberTitle', 'off');
   end;
   figure(hf);
   clf;
  
   for i = 1:length(sf)
      [tx, vx]=trace_ave2(sf(i), pflag);
      if(i == 1)
         v=zeros(length(sf), length(tx));
         t=zeros(length(sf), length(tx));
      end;
      v(i,:) = vx';
      t(i,:) = tx;
   end;
      [hp, err] = physplot(hf, 'A', 'Average', t, v, [0.1, 0.1, 0.6, 0.8], ...
         [0, max(max(t)), min(min(v)), max(max(v))], -60);
 	drawnow  
end;
%catch
%  watchoff;
%   QueMessage('Error in Current-clamp Analysis routine', 1);
%end;


function [ta, vmean] = trace_ave2(sf, plot_flag)

global VOLTAGE CURRENT DFILE
global CONTROL

err = 0;

QueMessage('Trace Average - Starting');
dat = [];
time = [];
ivresult = []; % initialize it
do_spike = 0; % assume we do not want it
[DFILE, err] = analysis_setup(DFILE, sf);

if(err ~= 0)
   return;
end;

% abstract general information
spike_thresh=number_arg(CONTROL(sf).thresh);
QueMessage(sprintf('CCIV: Spike Threshold: %7.2f mV', spike_thresh));
protocol=deblank(lower(CONTROL(sf).protocol));

[records,pts]=size(CURRENT);

[RL, err] = record_parse(CONTROL(sf).reclist);
frec=min(RL);
lrec=max(RL);

QueMessage('EPSP_TC - time base');
% compute the time base for plotting (time is [rec, npts] array for EACH record)
time=make_time(DFILE);
tmax=max(max(time));
k = find(DFILE.rate < 1);
DFILE.rate(k) = 50;
RATES = (DFILE.rate .* DFILE.nr_channel)/ 1000; % array of sampling rates, convert to msec


% always compute the time of each stimulus, in seconds from the start of the data
QueMessage('EPSP_TC - time base 2');
if(DFILE.mode >= 5)
   wz=DFILE.ztime;
   w=find(diff(wz) < -12*60*60); % correct for possibility that someone actually records the data across midnight rollover.... (yes, it happened. 5/16/01 with Huijie's data.)
   if(~isempty(w))
      wz(w+1:end)=wz(w+1:end)+24*60*60;
   end;
   zt = (wz-wz(1))/(60);
   cond_baseline = 5; % 5 min baseline
else
   zt = (DFILE.ztime-DFILE.ztime(1))/(60*1000);
   cond_baseline = 3;
end;
TM=zt;

QueMessage('EPSP_TC - time base 3');
if(DFILE.mode >= 5)
   wz=DFILE.ztime;
   w=find(diff(wz) < -12*60*60); % correct for possibility that someone actually records the data across midnight rollover.... (yes, it happened. 5/16/01 with Huijie's data.)
   if(~isempty(w))
      wz(w+1:end)=wz(w+1:end)+24*60*60;
   end;
   ZT = (wz-wz(1))/(60);
   cond_baseline = 5; % 5 min baseline
elseif length(DFILE.ztime) == length(RL)
   ZT = (DFILE.ztime-DFILE.ztime(1))/(60*1000);
   cond_baseline = 3;
else
   ZT = ones(length(RL), 1)*DFILE.cycle;
   cond_baseline = 5;
end;

% access ZT with the record number (1..n) to get the corresponding ztime
% have to get the records first

% Now,get the times when the valves switched (if any...)
% and generate periods with valve 1, 2, 3 or 4.
QueMessage('EPSP_TC - Valves...');
p=datac('getnote'); % read the current notefile information.
t_sw_valve=[];
n_valve=[];
TL=[];
VL=[];
if(~isempty(p) & length([p.proto]) > 0) % there should be some, but if not, don;'t do much
   % first set of arrays are immediate representations.
   sw_valve=[1 [p(find(diff([p.valve])~=0)+1).frec]];	% valve switch list (records)
   if(length(sw_valve) > 1)
      n_valve=[1 p(find(diff([p.valve])~=0)+1).valve]; % which valve...
      t_sw_valve=ZT(sw_valve); % don't forget offset from start of data..
      % now make long time arrays to match the other arrays.
      for i = frec:lrec
         TL(i)=ZT(i);
         for j=1:length(sw_valve)
            if(i >= sw_valve(j))
               VL(i)=n_valve(j);
            end;
         end;
      end;
   end;
end;

% put up a dialog to select the records in the traces
tmax = ceil(max(TM));
if(tmax < 12)
   part = 0.5;
else
   part = 0.75;
end;

dlgTitle = sprintf('Selection for Averaging: Time window is %7.2f to %7.2f', min(TM), max(TM));
prompt={'dt1','dt2', 'frequency', 'measurewin'};
def={sprintf('%7.2f', part*tmax), sprintf('%7.3f', tmax), '5', '55'};
dlgTitle='Selection parameters for STDP';
lineNo=1;
answer=inputdlg(prompt,dlgTitle,lineNo,def);
t1 = number_arg(answer{1});
t2 = number_arg(answer{2});
tt1 = number_arg(answer{3});
tt2 = number_arg(answer{4});


r1 = min(find(TM>= t1));
r2 = max(find(TM<= t2));
if(isempty(r1) | isempty(r2))
   QueMessage('Trace Average: window no good?');
   return;
end;

rate = mean(RATES);
tr1 = floor(tt1/rate);
tr2 = floor(tt2/rate);

vmean = mean(VOLTAGE(r1:r2, tr1:tr2))';
ta = time(1,tr1:tr2);
return;



